PROCEEDINGS OF THE 31"* ICRC, LODZ 2009 



1 



Ultra high energy cosmic rays: subluminal and superluminal 

shocks 

Athina Meli*, Julia K. Becker^^, John Quenby-f 



* Erlangen Center for Astroparticle Physics, Universitdt Erlangen-NUrnberg, Germany 
^ Institutionen for Fysik, Goteborgs Universitet, 41296 Goteborg, Sweden 
High Energy Physics Group, Blackett Laboratory, Imperial College London, Prince Consort Road SW7 2BZ, UK 



o 
o 

(N 
>>■ 



(N 

m 

X 
o 



> 

in 
o 
a^ 
o 



Abstract. Diffusive sliock acceleration is invoked 
to explain non-thermal particle acceleration in Su- 
pernova Remnants, Active Galactic Nuclei (AGN) 
Jets, Gamma ray Bursts (GRBs) and various large 
scale cosmic structures. The importance of achieving 
the highest observed particle energies by such a 
mechanism in a given astrophysical situation is a 
recurring theme. In this work, shock acceleration 
in relativistic shocks is discussed, mostly focusing 
on a numerical study concerning proton acceleration 
efficiency by subluminal and superluminal shocks, 
emphasising on the dependence of the scattering 
model, bulk Lorentz factor and the angle between 
the magnetic field and the shock flow. We developed 
a diffuse cosmic ray model based on the study 
of different shock boost factors, which shows that 
spectra from AGN fit current observations of ultra 
high energy cosmic rays, above 5.7 x 10^° GeV, 
much better than GRBs, indicating that AGN are 
the primary candidates to explain the UHECR flux. 
Recent Fermi observations of GRB090816c indicate 
very flat spectra which are expected within our model 
predictions and support evidence that GRB particle 
spectra can be flat, when the shock Lorentz factor is 
of order - 1000. 

Keywords: relativistic shocks, acceleration, AGN, 
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I. Introduction 

Active Galactic Nuclei (AGN) and Gamma Ray 
Bursts (GRBs) seem to be the two most promising 
source candidates for the production of charged ultra 
high energy cosmic rays (UHECRs). Work in the late 
1970s by a number of authors, e.g. H, lH, lfT4l . 
basing their ideas on the original Fermi acceleration 
mechanism, i.e. 121, IS], established the basic mecha- 
nism of particle diffusive acceleration in non-relativistic 
shocks. Since then, considerable analytical and numer- 
ical investigations have been performed, but questions 
remain concerning details of the acceleration mechanism 
at highly relativistic shock speeds and the consequent 
production and origin of UHECRs. 



II. The microphysics of shocks - numerical 

APPROACH 

In this two-fold work we firstly, present a series of 
Monte Carlo simulations with the aim to provide a 
more refined determination of the possible accelerated 
particle shock spectra that can result and secondly, we 
develop a diffuse cosmic ray (CR) model to explain the 
UHECR observed spectrum. A Monte Carlo numerical 
approach solves the well known Boltzmann transport 
equation which depends on the assumptions that the 
collisions represent diffusive scattering in pitch angle 
and that the scattering is elastic in the fluid frame where 
there is no residual electric field, taking into account 
the theoretical result that Alfven waves are limited to 
Va c/\/3- Then, the first order Fermi (diffusive) 
acceleration mechanism is simulated by following the 
test particles' guiding centres and allowing for numerous 
pitch angle scatterings in interaction with the assumed 
magnetised media, while at each shock crossing the 
particles gain an amount of energy determined by a 
Lorentz transformation of a reference frame. Standard 
theory assumes conservation of the first adiabatic in- 
variant in the so called de Hoffmann-Teller (HT) frame 
in order to determine reflection or transmission of the 
particles. Trajectory integration calculations, see ID, 
|5 |, giving the phase dependence of the probability of 
transmission as a function of phase and pitch angle, 
showed that the reflection percent plotted against pitch 
angle, never varied more than 20% from the mean value. 
In the relativistic shock situation, anisotropy renders the 
input to the shock from upstream, very anisotropic in 
pitch angle, but as discussed in |6|, it is an acceptable 
approximation to randomise phase before transforming 
to the HT frame and then to use the adiabatic invariant 
to decide on reflection or transmission. Therefore, for 
our simulations the initial injection energy (7) of the 
particles is taken to be 7 ~ (F + 10), where F is the 
shock Lorentz boost factor. 

For the oblique shock cases, especially studied here 
(for a detailed description of the shock kinematics see 
[;8l), we apply a pitch angle scatter [1/F < 56 < 
10/r,(/) S (0, 27r)] which approximately corresponds 
to a situation with a power spectrum of scattering 
waves, P{k)/ B'^ = 5/4^2 • F^^ . fc-i and with the 
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neglect of cross-field diffusion, right up to the shock 
interface. Provided the field directions encountered are 
reasonably isotropic in the shock frame, we know that 
tanT/^i = T^^ tanil^NSH ^ ~ V'l where '1' 

and 'NSH' refer to the upstream and normal shock 
frames respectively and ^ is the angle between the 
magnetic field and the shock normal. The concentration 
of field vectors close to the x-axis in the upstream 
fluid frame allows a reasonable probability of finding 
a HT frame with a boost along the negative y-axis less 
than c. Making this boost then yields an upstream HT 
frame inclination, ia,iiipHT,i — Tht,! tan?/'i. While 
all particles are allowed to cross from downstream to 
upstream, only particles with a critical HT frame pitch 
angle, 6^, given by 6 c = arcsin( )°'^ are allowed 
to cross upstream to downstream and conservation of 
the first adiabatic invariant is used to determine the new, 
downstream pitch angle |9|. The pitch angle is measured 
in the local fluid frame, while the value Xi gives the 
distance of the particles to the shock front, where the 
shock is assumed to be placed at a; = 0. A compression 
ratio of 3 is used and although some MHD conditions 
favour a value of 4, the study of |11| does not find a 
substantial difference between the simulation results for 
these two cases. Away from the shock, the guiding centre 
approximation is used so that a test particle moving a 
distance, d, along a field line at ijj to the shock normal, in 
the plasma frame has a probability of collision within d 
given by P{d) = 1— cxp(— d/ A) = R, where the random 
number i? is < i? < 1. Weighting the probability 
by the current in the field direction /i (i.e. cos 9) yields 
d = —XfilnR. 

For a detailed description on the numerical method see 
El. 

HI. Results - Discussion 

1 ) Superluminal shock spectra: We follow the helical 
trajectory of the particle until it intersects the shock 
front (applying pitch angle scatter [l/F < 56 < 
10/r, (j) e (0, 2tt)] as mentioned above). Simulation runs 
demonstrated that the results were almost independent 
of tjj. Shock frame spectra are illustrated in figure 1 
with a simulation run for an arbitrarily chosen, typically 
large inclination angle, ip — 76°, employing a range of 
shock boost factors. We see that superluminal relativistic 
shocks are not efficient accelerators for very high energy 
particles and are unlikely to contribute to observable 
effects as we will see lateron. These conclusions concur 
with the work of [12J. 

2) Subluminal shock spectra: Particle spectra pro- 
duced in relativistic subluminal shocks with the numeri- 
cal method of |8| have been calculated for three different 
inclination angles in the shock frame, — 23°, 33° and 
43°. Shown here (figure 2) are particle spectra averaged 
over the three inclination angles at a particular F (where 
F = 10, 100, 300, 500, 1000). The averaged values give a 
more realistic estimate of the diffuse particle flux (to be 
discussed later) from extragalactic sources, since a range 
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Fig. I. Superluminal, relativistic spectra at t/i = 76deg. Shock 
boost factors are varied between F = 10, 100, 300, 500, 1000. Spectra 
for different inclination angles ip are comparable. In the region of 
efficient acceleration, the spectra approximately follow power-laws 
with spectral indices lying between ~ 2.0 — 2.3. 
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Fig. 2. Subluminal spectra averaged over three shock angles (-0 = 
23°, 33° and 43°) for different T: F = 10,20,30 is displayed 
in the first row, in the middle, F = 100, 300, 500 is shown and 
F = 700, 900, 1000 is the bottom row. The black circles in each 
graph represent the simulation result. The red (solid) lines show the 
single power-law for comparison and the spectral indexes (a) values 
are shown accordingly. 



of angles is more likely to occur in the AGN and GRB 
shocks. A power-law fit to the simulated spectra between 
10^-^ GeV and lO^O-^ GeV is taken. We concentrate 
on the highest energies because observation of particle- 
induced air showers indicate that at energies between 
10^^ GeV and W^-^ GeV, the origin of the charged CR 
flux is extragalactic and that this flux is distinct from a 
dominant galaxy produced component at lower energies. 
At even higher energies, the spectra will be modified in 
practise by the absorption of protons due to interactions 
with the cosmic microwave background. 
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Fig. 3. The maximum predicted diffuse flux from GRBs and AGN 
with superluminal shock fronts (dashed line) and subluminal shocks 
(solid blue (upper) and black (lower) lines). For subluminal sources, 
the upper line is a pure AGN-produced spectrum, the lower line 
represents a pure GRB spectrum. The flux is compared to the measured 
CR spectrum. In the case of superluminal sources, 50% is assumed to 
come each from GRB and 50% from AGN. 
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Fig. 4. Spectrum of UHECRs multiplied by E^'^. Data points 
from Auger, 1 16| and HiRes, 1 19|. The solid lines represent the same 
predictions as presented in Fig. [5] Auger data have been renormalised 
at 10^" '^ GeV to HiRes data and the calculated spectra have also been 
normalised to HiRes data. The data can be described well by a pure 
AGN spectrum (blue (lower) hne) within experimental uncertainties 
and test particle acceleration accuracy. The red (middle) line is a 
mixture of 50% GRBs contribution and 50% AGN contents, the black 
(upper) line is a pure GRB spectrum. 



Regarding relativistic shocks and their characteristics, 
ll2n investigating parallel shocks up to F = 30 found 
an increase in structure in the spectral shape and a 
decrease in slope as V increased with a dependence 
^-1.26 at r = 30. These trends are found in our 
work, extending to far higher T factors and for a 
more general set of subluminal shock inclination angles. 
Furthermore, [22J employed pitch angle scattering and 
varying cross-field diffusion coefficients and found that 
at low r, steep spectra occurred at large inclination 
angles but all values of these parameters seemed to 
produce spectral slopes of -2.2 at F = 243. In contrast 
ifTOll found spectra flatter than E^'^ for parallel shocks 



as F ^ 1000. In addition, related work by |l23l wifli 
wave spectra P{k) ^ found spectra flatter 

than E^"^ with noticeable spectral structure in inclined, 
subluminal shocks at upstream velocities of 0.5c and a 
weakly perturbed field. Their trajectory integrations took 
into account cross field diffusion. There is an agreement 
that at very high inclinations, significant acceleration 
above that due to a single shock cycle is ruled out. 
One realises that there seems to be a developing consen- 
sus that high F subluminal shocks result in flatter spectra 
than E~^. In il5i . it is shown that GRBs can have 
flat relativistic electron spectra. A number of varying 
power-law slopes are also consistent with radio data on 
the electron spectra injected at terminal hotspots in the 
lobes of powerful FR-II radio galaxies i.e. |18| among 
others. Recenlty, Fermi observations of GRB090816c 
||29l indicate very flat spectra, to be discussed lateron. 

3) The diffuse spectrum of CRs: The source proton 
spectra derived above can be translated into an expected 
diffuse proton flux from astrophysical sources by folding 
the spectra with the spatial distribution of the sources. 
In our calculations, adiabatic energy losses are taken 
into account and it is also assumed that both AGN and 
GRBs follow the star formation rate to determine the 
number density evolution with comoving volume, see 
for example, (20\ in the case of AGN and [24| in the 
case of GRBs. Also the absorption of protons due to the 
interactions with the CMBR at Ep > 5- 10^^ eV, as was 
recently confirmed by the Auger experiment, see ifTSl . is 
taken into account. Furthermore, because the calculated 
particle spectra are given in arbitrary units, normalisation 
of the overall spectrum as measured at Earth is achieved 
using observation. 

(i) In the case of superluminal sources, normalisation 
of the expected signal follows from the most restrictive 
upper limit on the neutrino signal from extraterrestrial 
sources given by the AMANDA experiment ( 1261 ). 
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f dE,y ^ 1 .1 - s sr cm^ ■ 
(ii) In the case of subluminal sources, using neutrino flux 
limits leads to an excess above the observed spectrum of 
charged CRs, since the limits are not stringent enough 
yet. Instead, the measured CR spectrum above the 'an- 
kle' is used to estimate the contribution from subluminal 
sources. The CR energy flux above the ankle is given 
by e.g. L30J. 

It is expected that the contribution comes from a com- 
bined signal from AGN and GRBs. It is assumed that the 
fraction of UHECRs coming from AGN, contributes a 
fraction < a; < 1. Therefore, the fraction of UHECRs 
from GRBs is (1 — x). Thus, the total spectrum as 
observed at Earth is given as: 

^ = ^(^p(-)) + (i - ■ 

The minimum redshift is set to z = 0.001. This 
excludes only the closest AGNs which contribute less 
than 1% according to newest results of ||251 . Il27l . The 
maximum redshift is taken to be z^ax — 7. As the main 
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contribution comes from redshifts of z ~ 1 — 2 due to 
the high number of sources at these redshifts, the exact 
values of the integration Hmits are not crucial. 
The diffuse spectrum from superluminal and subluminal 
shock sources (dashed, solid lines) as could be measured 
at the Earth is shown in figure 3 and the diffuse spectrum 
multiplied by E'^-'' , focusing on energies > 10^ GeV for 
subluminal shock sources, is shown in figure 4. One sees 
a very good model agreement (blue (lower) line of figure 
4). 

As aforementioned, Fermi observed the GRB090816c 
with both instruments, GBM and LAT ||29l . This means 
that the spectral behaviour of this burst is studied 
between 10 keV and 10 GeV, i.e. over 6 orders of 
magnitude. The observed spectrum is believed to arise 
from electron synchrotron radiation. The observation of 
the highest photons gives a limit on the boost factor, 
as no cutoff due to pair production was observed. The 
minimum boost factor is varying for the different time 
bins and it ranges between Train ^ 600 — 1100. The 
photon spectrum at energies below the break energy, 
where it is expected to be unaffected by absorption 
effects, shows a spectral behavior between v^^-^ — v^^ . 
In the fast cooling regime, this implies primary electron 
spectra of E^^-^ — E^"^-^ . These very flat spectra are 
expected within our model and supporting evidence 
that GRB particle spectra can be flat. Note however, 
that electrons easily interact and primary spectra are 
therefore easily modified. Hence, it is not clear yet 
if the observed spectrum is indeed unmodified. The 
observation of neutrino spectra from GRBs can help to 
solve that matter, as neutrinos freely escape from the 
source, see e.g. Il28l . 

IV. Summary & conclusions 

We have presented Monte Carlo simulation studies 
of the acceleration of test particles in relativistic, 
subluminal and superluminal shock environments (i.e. 
AGN, GRBs) and the resulting CR spectra were used 
to calculate a diffuse CR contribution. Our results can 
be summarised as follows: 

1. Subluminal shock acceleration was studied with a 
pitch scattering angle. The resulting spectral slopes were 
roughly independent of inclination angle, though some 
details of the features were different. A dependence 
of the spectral index ap on the shock boost factor 
r was found, leading to spectra of ap ^ 2.0 — 2.1 
for mildly-relativistic shocks of F ^ 10 — 30, but 
producing much harder spectra (1.0 < ap < 1.5) for 
highly-relativistic shocks, 100 < F < 1000. The particle 
spectra arising from relativistic shocks in GRB with 
very high boost factors between 100 < F < 1000, 
have spectra flatter than iJp^^'^. The above findings 
are supported by the work at lower F of ll23l and | |2TT |. 
Observational evidence, see [15], regarding irregular 
and flat spectra from GRBs may be explained by the 
spectra we present. 

2. Superluminal shocks, applying the same pitch angle 



scatter as above, seem only efficient in accelerating 
CRs up to Ep ^ 10^ GeV, resulting in spectral indices 
of Op ^ 2.0 — 2.3. On the other hand, subluminal 
shocks are more efficient and able to accelerate CRs 
up to Ep - 10^2 GeV, factors of 10^^" above the 
particle injection energy. Extragalactic, superluminal 
shocks energy density at low energies may be quite 
high compared to the observed flux of UHECRs, but 
hidden by galactic cosmic rays. Nevertheless they could 
be good candidates for the production of high energy 
neutrinos and photons, see 

3. Diffuse UHECR flux contributions of AGN and 
GRBs was discussed. For superluminal sources, such 
contributions can be excluded by using current neutrino 
flux limits to normalise the spectrum. In the case of 
subluminal sources, the spectrum is normalised to 
the CR flux above the ankle, E^in = 10^^ GeV. By 
applying AGN (with shock boost factors of F = 10) as 
potential UHECRs sources, the spectrum fits very well 
the current CR data within experimental uncertainties. 

Recent Fermi observations of GRB090816c indicate 
very flat spectra which are expected within our model 
predictions and support evidence that GRB particle 
spectra can be flat, when the shock Lorentz factor is 
of order - 1000. 
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